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This chapter introduces the fracture nucleation process, their (extreme) statistics in disordered 
solids, in fiber bundle models, and in the two fractal overlap models of earthquake. 



I. INTRODUCTION 
Models of Fracture in Disordered Solids and Statistics 



If one applies tensile stress on a solid, the solid elongates and gets strained. The stress (a) - strain (e) relation 
is linear for small stresses (Hooke's law) after which nonlinearity appears, in most cases. Finally at a critical stress 
o f, depending on the material, amount of disorder and the specimen size etc., the solid breaks into pieces; fracture 
occurs. In the case of brittle solids, the fracture occurs immediately after the Hookean linear region, and consequently 
the linear elastic theory can be applied to study the essentially nonlinear and irreversible static fracture properties of 
brittle solids p|. With extreme perturbation, therefore, the mechanical or electrical properties of solids tend to get 
destabilised and failure or breakdown occurs. In fact, these instabilities in the solids often nucleate around disorder, 
which then plays a major role in the breakdown properties of the solids. The growth of these nucleating centres, in 
turn, depends on various statistical properties of the disorder, namely the scaling properties of percolating structures, 
its fractal dimensions, etc. These statistical properties of disorder induce some scaling behaviour for the breakdown 
of the disordered solids 0, Q . 

Obviously with more and more random voids, the linear response of e.g. the modulus of elasticity Y (say, the 
Young's modulus) of the solid decreases. So also does the breaking strength of the material: the fracture strength 
1 Of of the specimen. For studying most of these mechanical (elastic) breakdown problems of randomly disordered 
\& , solids, one can take the lattice model of disordered solids. In these lattice models, a fraction p of the bonds (or sites) 
0^ ' are intact, with the rest (1 — p) being randomly broken or cut. Fluctuations in the random distribution give rise to 
random clusters of springs inside the bulk, for which the statistics is well developed 4], and one can investigate the 
effect of the voids or impurity clusters on the ultimate strength of the (percolating) elastic network of the bulk solid 
0, Q • One can also consider and compare the results for failure strength of the solids with random bond strength 
distribution, as for example, in random fiber bundle models (see e.g. 00,0)- 

As is well known, the initial variations (decreases) of the linear responses like the elastic constant Y so the breakdown 
strengths aj is analytic with the impurity (dilution) concentration. Near the percolation threshold Q p c , up to (and 
at) which the solid network is marginally connected through the nearest neighbour occupied bonds or sites and below 
which the macroscopic connection ceases, the variations in these quantities with p are expected to become singular; the 
leading singularities being expressed by the respective critical exponents. The exponents for the modulus of elasticity 
Y ~ Ap Tc Ap = (p — p c )/Pc are well-known and depend essentially on the dimension d of the system (see e.g., £j|). 
J> . One kind of investigation searches for the corresponding singularities for the essentially nonlinear and irreversible 
properties of such mechanical breakdown strengths for p near the percolation threshold p c : to find the exponent Tf 
for the average fracture stress 

tjf ~ Ap Tf for p near p c . Very often, one maps [j| the problem of breakdown to the corresponding linear problem 
(assuming brittleness up to the breaking point) and then derives the scaling relations giving the breakdown 
exponents Tf in terms of the linear response exponent (T e ) and other lattice statistical exponents (see next section). 

Unlike that for the 'classical' linear responses of such solids, the extreme nature of the breakdown statistics, 
nucleating from the weakest point of the sample, gives rise to a non-self-averaging property. We will discuss (in the 
next section) these distribution functions F(a), giving the cumulative probability of failure of a disordered sample of 
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linear size L. We show that the generic form of the function F{a) can be either the Weibull [2] form 



or the Gumbel [2] form 



F(a) ~ 1 — exp 



F(a) ~ 1 — exp 
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where A(p) is determined by the linear response like the elasticity of the disordered solid by some other lattice 
statistical quantity etc. and (f> is an exponent discussed in the next section. 

In another kind of model for disordered systems, a loaded bundle of fibers represents the various aspects of fracture 
process through its self-organised dynamics. The fiber bundle model study was initiated by Pierce j6j in the context 
of testing the strength of cotton yarns. Since then, this model has been studied from various points of view. Fiber 
bundles are of two classes with respect to the time dependence of fiber strength: The 'static' bundles contain fibers 
whose strengths are independent of time, whereas the 'dynamic' bundles are assumed to have time dependent elements 
to capture the creep rupture and fatigue behaviors. For simplicity, we will discuss here the 'static' fiber bundle models 
only. According to the load sharing rule, fiber bundles are being classified into two groups: Equal load-sharing (ELS) 
bundles or democratic bundles and local load-sharing (LLS) bundles. In democratic or ELS bundles, intact fibers bear 
the applied load equally and in local load-sharing bundles the terminal load of the failed fiber is given equally to all the 
intact neighbors. The classic work of Daniels on the strength of the static fiber bundles under equal load sharing 
(ELS) assumption initiated the probabilistic analysis of the model (see e.g., @). The distribution of burst avalanches 
during fracture process is a marked feature of the fracture dynamics and can be observed in ultrasonic emissions 
during the fracture process. It helps characterizing different physical systems along with the possibility to predict the 
large avalanches. From a nontrivial probabilistic analysis, one gets power law distribution of avalanches for static 
ELS bundles, whereas the power law exponent observed numerically for static LLS bundles differs significantly. This 
observation induces the possibility of presenting loaded fiber bundles as earthquake models (see Sec. 3). The phase 
transition [8| and dynamic critical behavior of the fracture process in such bundles has been established through 
recursive formulation 0, of the failure dynamics. The exact solutions [lj} of the recursion relations suggest 
universal values of the exponents involved. Attempt has also been made ^lj to study the ELS ans LLS bundles from 
a single framework introducing a 'range of interaction' parameter which determines the load transfer rule. 



B. Earthquake Models and Statistics 



The earth's solid outer crust, of about 20 kilometers in average thickness, rests on the tectonic shells. Due to 
the high temperature-pressure phase changes and the consequent powerful convective flow in the earth's mantle, at 
several hundreds of kilometers of depth, the tectonic shell, divided into a small number (about ten) of mobile plates, 
has relative velocities of the order of a few centemeters per year [l2l[T^ |. Over several tens of years, enormous elastic 
strains develop sometimes on the earth's crust when sticking (due to the solid-solid friction) to the moving tectonic 
plate. When slips occur between the crust and the tectonic plate, these stored elastic energies are released in 'bursts', 
causing the damages during the earthquakes. Because of the uniform motion of the tectonic plates, the elastic strain 
energy stored in a portion of the crust (block), moving with the plate relative to a 'stationary' neighbouring portion 
of the crust, can vary only due to the random strength of the solid-solid friction between the crust and the plate. The 
slip occurs when the accumulated stress exceeds the frictional force. 

As in fracture (in fiber bundle model in particular), the observed distribution of the elastic energy release in 
various earthquakes seems to follow a power law. The number of earthquakes N(m), having magnitude in the Richter 
scale greater than or equal to to, is phenomenologically observed to decrease with m exponentially. This gives the 
Gutenberg- Richter law [lj 

\nN{m) = constant — am, 

where a is a constant. It appears that the amount of energy e released in an earthquake of magnitude to is 

related to it exponentially: 

lne = constant + b to, 

where b is another constant. Combining therefore we get the power law giving the number of earthquakes N(e) 
releasing energy equal to e as 



N(e) ~ e- 



(3) 
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with a — a/b. The observed value of the exponent (power) a in is around unity (see e.g., 

Several laboratory and computer simulation models have recently been proposed Q to capture essentially the 
above power law in the earthquake energy release statistics. In a very successful table-top laboratory simulation 
model of earthquakes, Burridge and Knopoff [3] took a chain of wooden blocks connected by identical springs to the 
neighbouring blocks. The entire chain was placed on a rigid horizontal table with a rough surface, and one of the end 
blocks was pulled very slowly and uniformly using a driving motor. The strains of the springs increase due to the 
creep motions of the blocks until one or a few of the blocks slip. The drops in the elastic energy of the chain during 
slips could be measured from the extensions or compressions of all the springs, and could be taken as the released 
energies in the earthquake. For some typical roughness of the surfaces (of the blocks and of the table), the distribution 
of these drops in the elastic energy due to slips indeed shows a power law behaviour with a ~ 1 in |J3J. 

A computer simulation version of this model by Carlson and Langer [l5l | considers harmonic springs connecting 
equal mass blocks which are also individually connected to a rigid frame (to simulate other neighbouring portions of 
the earth's crust not on the same tectonic plate) by harmonic springs. The entire system moves on a rough surface 
with nonlinear velocity dependent force (decreasing to zero for large relative velocities) in the direction opposite to 
the relative motion between the block and the surface. In the computer simulation of this model it is seen that the 
distribution of the elastic energy release in such a system can indeed be given by a power law like @, provided the 
nonlinearity of the friction force, responsible for the self-organisation, is carefully chosen |15| . 

The lattice automata model of Bak et al represent the stress on each block by a height variable at each lattice 
site. The site topples (the block slips) if the height (or stress) at that site exceeds a preassigned threshold value, and 
the height becomes zero there and the neighbours share the stress by increasing their heights by one unit. With this 
dynamics for the system, if any of the neighbouring sites of the toppled one was already at the threshold height, the 
avalanche continues. The boundary sites are considered to be all absorbing. With random addition of heights at a 
constant rate (increasing stress at a constant rate due to tectonic motion), such a system reaches its self-organised 
critical point where the avalanche size distributions follow a natural power law corresponding to this self-tuned critical 
state. Bak et al identify this self-organised critical state to be responsible for the Gutenberg-Richter type power 
law. All these models are successful in capturing the Gutenberg-Richter power law, and the real reason for the self- 
similarity inducing the power law is essentially the same in all these different models: emergence of the self-oganised 
critical state for wide yet suitably chosen variety of nonlinear many-body coupled dynamics. In this sense all these 
models incoporate the well-established fact of the stick-slip frictional instabilities between the earth's crust and the 
tectonic plate. It is quite difficult to check at this stage any further details and predictions of these models. 

While the motion of the tectonic plate is surely an observed fact, and this stick-slip process should be a major 
ingredient of any bonafide model of earthquake, another established fact regarding the fault geometries of the earth's 
crust is the fractal nature of the roughness of the surfaces of the earth's crust and the tectonic plate. This latter 
feature is missing in any of these models discussed above. In fact, the surfaces involved in the process are results of 
large scale fracture seperating the crust from the moving tectonic plate. Any such crack surface is observed to be 
a self-similar fractal, having the self-affine scaling property z(Xx, Ay) ~ X^z(x, y) for the surface coordinate z in the 
direction perpendicular to the crack surface in the (x,y) plane 2\. Various fractographic investigations indicate a 
fairly robust universal behaviour for such surfaces and the roughness exponent £ is observed to have a value around 
0.80 — 0.85 (with a possible crossover to £ ~ 0.4 for slow propagation of the crack-tip) 0, Q- This widely observed 
scaling property of the fracture surfaces also suggests that the fault surfaces of the earth's crust or the tectonic plate 
should have similar fractal properties. In fact, some investigators of the earthquake dynamics have already pointed 
out that the fracture mechanics of the stressed crust of the earth forms self-similar fault patterns, with well-defined 
fractal dimensionalities near the contact areas with the major plates Based on these observations regarding the 
earthquake faults, we have developed a 'two-fractal overlap' model [13, 13 discussed later. 



II. FRACTURE STATISTICS OF DISORDERED SOLIDS 
A. Griffith Energy Balance and Brittle Fracture Strength of Solids 

One can easily show (see e.g. 0,0) that in a stressed solid, the local stress at sharp notches or corners of the 
microcrack can rise to a level several times that of the applied stress. This indicetes how the microscopic cracks 
or flaws within a solid might become potential sources of weakness of the solid. Although this stress concentration 
indicates clearly where the instabilities should occur, it is not sufficient to tell us when the instability does occur and 
the fracture propagation starts. This requires a detailed energy balance consideration. 

Griffith in 1920, equating the released elastic energy (in an elastic contin- uum) with the energy of the surface newly 
created (as the crack grows) , arrived at a quantitative criterion for the equilibrium extension of the microcrack already 
present within the stressed material [19J . We give below an analysis which is valid effectively for two-dimensional 
stressed solids with a single pre-existing crack, as for example the case of a large plate with a small thickness. Extension 
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FIG. 1: A portion of a plate (of thickness w) under tensile stress a (Model I loading) containing a linear crack of length 21. For 
a further growth of the crack length by 2dl, the elastic energy released from the annular region must be sufficient to provide 
the surface energy 4Fu;d/ (extra elastic energy must be released for finite velocity of crack propagation). 

to three-dimensional solids is straightforward. 

Let us assume a thin linear crack of length 21 in an infinite elastic continuum subjected to uniform tensile stress 
I perpendicular to the length of the crack (see Fig. I). Stress parallel to the crack does not affect the stability of 
the crack and has not, therefore, been considered. Because of the crack (which can not support any stress field, at 
least on its surfaces), the strain energy density of the stress field (a 2 /2Y) is perturbed in a region around the crack, 
having dimension of the length of the crack. We assume here this perturbed or stress-released region to have a circular 
cross-section with the crack length as the diameter. The exact geometry of this perturbed region is not important 
here, and it determines only an (unimportant) numerical factor in the Griffith formula (see e.g. pj). Assuming 
therefore half of the stress energy of the annular or cylindrical volume, having the internal radius / and outer radius 
I + dl and length w (perpendicular to the plane of the stress; here the width w of the plate is very small compared to 
the other dimensions), to be released as the crack propagates by a length dl, one requires this released strain energy 
to be sufficient for providing the surface energy of the four new surfaces produced. This suggests 



Here Y represents the Young's modulus of the solid and T represents the surface energy density of the solid, measured 
by the extra energy required to create unit surface area within the bulk of the solid. 

We have assumed here, on average, half of the strain energy of the cylindrical region having a circular cross-section 
with diameter 21 to be released. If this fraction is different or the cross-section is different, it will change only some 
of the numerical factors, in which we are not very much interested here. Also, we assume here linear elasticity up to 
the breaking point, as in the case of brittle materials. The equality holds when energy dissipation, as in the case of 
plastic deformation or for the propagation dynamics of the crack, does not occur. One then gets 



for the critical stress at and above which the crack of length 21 starts propagating and a macroscopic fracture occurs. 
Here A is called the stress-intensity factor or the fracture toughness. In fact, one can alternatively view the fracture 
occurring when the stress at the crack-tip (given by the stress intensity factor in exceeds the elastic stress limit 
for the medium. 

In a three-dimensional solid containing a single elliptic disk-shaped pla- nar crack perpendicular to the applied tensile 
stress direction, a straight- forward extension of the above analysis suggests that the maximum stress concentration 





(4) 
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would occur at the two tips (at the two ends of the major axis) of the ellipse. The Griffith stress for the brittle fracture 
of the solid would therefore be determined by the same formula Q , with the crack length 21 replaced by the length of 
the major axis of the elliptic planar crack. Generally, for any dimension therefore, if a crack of length I already exists 
in an infinite elastic continuum, subject to uniform tensile stress a perpendicular to the length of the crack, then for 
the onset of brittle fracture , Griffith equates (the differentials of) the elastic energy Ei with the surface energy E s : 



where Y represents the elastic modulus appropriate for the strain, T the surface energy density and d the dimension. 
Equality holds when no energy dissipation (due to plasticity or crack propagation) occurs and one gets 

^-AjA-n/FI (6) 

for the breakdown stress at (and above) which the existing crack of length I starts propagating and a macroscopic 
fracture occurs. It may also be noted that the above formula is valid in all dimensions (d > 2). 

For disordered solids, let us model the solid by a percolating system. As mentioned earlier, for the occupied 
bond/site concentration p > p c , the percolation threshold, the typical pre-existing cracks in the solid will have the 
dimension (I) of correlation length £ ~ Ap~ u and the elastic strength Y ~ Ap Ta Q. Assuming that the surface 
energy density T scales as with the backbone (fractal) dimension ds 0, equating E\ and E s as in ©, one gets 

(|l)f d ~£ dB - This gives 



with 



a f ~ (Apff 

T f = ±[T e + (d-d B )v] (7) 

for the 'average' fracture strength of a disordered solid (of fixed value) as one approaches the percolation threshold. 
Careful extensions of such scaling relations Q and rigorous bounds for Tt has been obtained and compared extensively 
in Hi- 



B. Extreme Statistics of The Fracture Stress 



The fracture strength af of a disordered solid does not have self-averaging statistics; most probable and the average 
o f may not match because of the extreme nature of the statistics. This is because, the 'weakest point' of a solid 
determines the strength of the entire solid, not the average weak points! As we have modelled here, the statistics of 
clusters of defects are governed by the random percolation processes. We have also discussed, how the linear responses, 
like the elastic moduli of such random networks, can be obtained from the averages over the statistics of such clusters. 
This was possible because of the self-averaging property of such linear responses. This is because the elasticity of a 
random network is determined by all the 'parallel' connected material portions or paths, contributing their share in 
the net elasticity of the sample. However, the fracture or breakdown property of a disordered solid is determined by 
only the weakest (often the longest) defect cluster or crack in the entire solid. Except for some indirect effects, most 
of the weaker or smaller defects or cracks in the solid do not determine the breakdown strength of the sample. The 
fracture or breakdown statistics of a solid sample is therefore determined essentially by the extreme statistics of the 
most dangerous or weakest (largest) defect cluster or crack within the sample volume. 

We discuss now more formally the origin of this extreme statistics. Let us consider a solid of linear size L, containing 
n cracks within its volume. We assume that each of these cracks have a failure probability fi(cr),i — 1,2, ... , n to 
fail or break (independently) under an applied stress a on the solid, and that the perturbed or stress-released regions 
of each of these cracks are seperate and do not overlap. If we denote the cumulative failure probability of the entire 
sample, under stress a, by F(a) then Q 



1 - F(<r) = jj(l - /*(<0) - ex P 



exp [-L d g{a)\ (8) 



where g(u) denotes the density of cracks within the sample volume L d (coming from the sum over the entire 
volume), which starts propagating at and above the stress level a. The above equation comes from the fact that the 
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FIG. 2: Schematic variation of failure probability F(a) with stress a for a disordered solid with volume Lf or (Li > Li). 

sample survives if each of the cracks within the volume survives. This is the essential origin of the above extreme 
statistical nature of the failure probability F(a) of the sample. 

Noting that the pair correlation g(l) of two occupied sites at distance I on a percolation cluster decays as 
exp (— l/£,(p)), and connecting the stress a with the length I by using Griffith's law that a ~ A, one gets 

g(a) ~ exp ^— ^373 ) for p — > p c . This, put in eqn. JHJ) gives the Gumbel distribution J2J) given earlier 0. If, 
on the other hand, one assumes a power law decay of g(l): g(l) ~ Z - ™, then using the Griffith's law QJ, one gets 
g(a) ~ (t)"\ giving the Weibull distribution (Q), from eqn. ©, where m — w/4> gives the Weibull modulus Q The 
variation of -F(c) with a in both the cases have the generic form shown in Fig. 2. F(a) is non-zero for any stress 
a > and its value (at any a) is higher for larger volume (L d ). This is because, the possibility of a larger defect 
(due to fluctuation) is higher in a larger volume and consequently, its failure probability is higher. Assuming F(af) 
is finite for failure, the most probable failure stress a/ becomes a decreasing function of volume if extreme statistics 
at work. 

The precise ranges of the validity of the Weibull or Gumbel distributions for the breakdown strength of disordered 
solids are not well established yet. However, analysis of the results of detailed experimental and numerical studies 
of breakdown in disordered solids seem to suggest that the fluctuations of the extreme statistics dominate for small 
disorder 0. Very near to the percolation point, the percolation statistics takes over and the statistics become self- 
averaging. One can argue |l9|. that arbitrarily close to the percolation threshold, the fluctuations of the extreme 
statistics will probably get suppressed and the percolation statistics should take over and the most probable breaking 
stress becomes independent of the sample volume (its variation with disorder being determined, as in Eqn.0, by 
an appropriate breakdown exponent). This is because the appropriate competing length scales for the two kinds 
of statistics are the Lifshitz scale In L (coming from the finiteness of the volume integral of the defect probability: 
L d (l — p) 1 finite, giving the typical defect size I ~ InL) and the percolation correlation length £. When £ < InL, 
the above scenario of extreme statistics should be observed. For £ > InL, the percolation statistics is expected to 
dominate. 

C. Failure Statistics in Fiber Bundles 

The fiber bundle (see Fig. consists of N fibers or Hook springs, each having identical spring constant k. The 
bundle supports a load W = Na and the breaking threshold (cr t / l ) i of the fibers are assumed to be different for different 
fiber (i). For the equal load sharing model we consider here, the lower platform is absolutely rigid, and therefore no 
local deformation and hence no stress concentration occurs anywhere around the failed fibers. This ensures equal load 
sharing, i.e., the intact fibers share the applied load W equally and the load per fiber increases as more and more 
fibers fail. The strength of each of the fiber (ath)i in the bundle is given by the stress value it can bear, and beyond 
which it fails. The strength of the fibers are taken from a randomly distributed normalised density p(cr t h) within the 
interval and 1 such that 

/ p(<Jth)d<Jth = 1. 
Jo 

The equal load sharing assumption neglects 'local' fluctuations in stress (and its redistribution) and renders the model 
as a mean-field one. 
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FIG. 3: The fiber bundle consists initially of N fibers attached in parallel to a fixed and rigid plate at the top and a downwardly 
movable platform from which a load W is suspended at the bottom. In the equal load sharing model considered here, the 
platform is absolutely rigid and the load W is consequently shared equally by all the intact fibers. 

i i 




FIG. 4: The simple model considered here assumes uniform density p(ath) of the fiber strength distribution up to a cutoff 
strength (normalized to unity). At any load per fiber level at at time t, the fraction a t fails and 1 — at survives. 

The breaking dynamics starts when an initial stress a (load per fiber) is applied on the bundle. The fibers having 
strength less than a fail instantly. Due to this rupture, total number of intact fibers decreases and rest of the (intact) 
fibers have to bear the applied load on the bundle. Hence effective stress on the fibers increases and this compels 
some more fibers to break. These two sequential operations, namely the stress redistribution and further breaking of 
fibers continue till an equilibrium is reached, where either the surviving fibers are strong enough to bear the applied 
load on the bundle or all fibers fail. 

This breaking dynamics can be represented by recursion relations in discrete time steps. For this, let us consider 
a very simple model of fiber bundles where the fibers (having the same spring constant k) have a white or uniform 
strength distribution p{uth) upto a cutoff strength normalized to unity, as shown in Fig. Ill CI p{&th) = 1 for 
< <Jth < 1 and = p{uth) — for a > tj t h- Let us also define Ut{<j) to be the fraction of fibers in the bundle that 
survive after (discrete) time step t, counted from the time t = when the load is put (time step indicates the number 
of stress redistributions). As such, U t {<r = 0) = 1 for all t and U t (a) = 1 for t = for any a; U t {<r) = U*(<j) ^ for 
t — * oo and a < a c , the critical or failure strength of the bundle, and Ut{cr) = for t — > oo if a > a c . 

Therefore U t (cr) follows a simple recursion relation (see Fig. Ill CJ|) 
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or, U t+1 = l-^-. (9) 

At the equilibrium state {Ut+i = U = U*), the above relation takes a quadratic form of U* : 

U* 2 - U* + a = 0. 

The solution is 

U*{a) = \±{af-af^;a f = \. 

Here cr/ is the critical value of initial applied stress beyond which the bundle fails completely. The solution with (+) 
sign is the stable one, whereas the one with (— ) sign gives unstable solution [10(. The quantity U*(a) must be real 
valued as it has a physical meaning: it is the fraction of the original bundle that remains intact under a fixed applied 
stress a when the applied stress lies in the range < a < af. Clearly, U*(0) = 1. Therefore the stable solution can 
be written as 

U*(a) = U*(a f ) + (a f - a) 1 / 2 ; U*(a f ) = \ and a f = ~. (10) 

For a > cr/ we can not get a real- valued fixed point as the dynamics never stops until Ut — when the bundle breaks 
completely. 

(a) At a < a/ 



It may be noted that the quantity U*{a) — U*(af) behaves like an ord er p arameter that determines a transition 

= U*(a)-U*(a f ) = (a f -af;P = ^-. (11) 



ierp 

from a state of partial failure (a < af) to a state of total failure (a > af) |10| : 

1 



2 

To study the dynamics away from criticality (a — > af from below), we replace the recursion relation by a 
differential equation 

dU U 2 -U + a 



dt U 

Close to the fixed point we write U t (a) = U*(a) + e (where e — > 0). This, following Eq. (10), gives [l"rT | 

e = U t (a)-U*{*)Kexp(-t/T), (12) 
where r = \ — cr) -1 / 2 + l] . Near the critical point we can write 

TK(a f -a)- a ;a=± (13) 

Therefore the relaxation time diverges following a power-law as a — > cr/ from below |Tfij . 

One can also consider the breakdown susceptibility \, defined as the change of U*(a) due to an infinitesimal 
increment of the applied stress a 



X = 



dU*(a) 



da 



from equation (10). Hence the susceptibility diverges as the applied stress a approaches the critical value cry = j. 
Such a divergence in x had already been observed in the numerical studies. 



(b) At a = aj 

At the critical point (cr = cr/), we observe a different dynamic critical behavior in the relaxation of the failure process. 
From the recursion relation JjJJl, it can be shown that decay of the fraction Ut(af) of unbroken fibers that remain 
intact at time t follows a simple power-law decay [Tol | : 

U t = \(l + ^ (15) 

starting from Uq = 1. For large t (i — > oo), this reduces to U t — 1/2 oc t~ 5 ; S = 1; a strict power law which is a robust 
characterization of the critical state. 
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FIG. 5: The fiber breaking strength distribution p(ath) 
behavior of the equal load sharing model. 



considered for studying elastic-plastic type nonlinear deformation 



1. Universality Class of The Model 

The universality class of the model has been checked ^} taking two other types of fiber strength distributions: (1) 
linearly increasing density distribution and (II) linearly decreasing density distribution within the (a t h) limit and 
1. One can show that while a j changes with different strength distributions (o/ = for case (I) and 0/ = 4/27 

for case II), the critical behavior remains unchanged: a = 1/2 = (3 = 7, 6 = 1 for all these equal load sharing models. 



2. Nonlinear Stress-Strain Relation for The Bundle 



One can now consider a slightly modified strength distribution of the equal load sharing fiber bundle, showing typical 
nonlinear deformation characteristics [71 llCj. For this, we consider an uniform density distribution of fiber strength, 
having a lower cutoff. Until failure of any of the fibers (due to this lower cutoff), the bundle shows linear elastic 
behavior. As soon as the fibers start failing, the stress-strain relationship becomes nonlinear. The dynamic critical 
behavior remains essentially the same and the static (fixed point) behavior shows elastic-plastic like deformation 
before rupture of the bundle. 

Here the fibers are elastic in nature having identical force constant k and the random fiber strengths distributed 
uniformly in the interval [ol, 1] with or, > 0; the normalised distribution of the threshold stress of the fibers thus has 
the form (see Fig. 5): 



P(°th) 



0. 

1 

1-0L 



0<(?th< &L 

, a L < a t h < 1 



(16) 



For an applied stress o < ol none of the fibers break, though they are elongated by an amount <5 = a j k. The 
dynamics of breaking starts when applied stress o becomes greater than ar,. Now, for a > ol the fraction of unbroken 
fibers follows a recursion relation (for p(a t h) as in Fig. 5): 



1 



F 
NUt 



1 



1 



1-OL 1-OL 



1 - — 



which has stable fixed points: 



U*(a) 



1 



2(1 -a L ) 



o 

1 

°7 



1/2 



1 



;°7 



4(1 - o L ) 



(17) 



(18) 



The model now has a critical point o^ = 1/[4(1 — ol)] beyond which total failure of the bundle takes place. The 
above equation also requires that ol < 1/2 (to keep the fraction U* < 1). As one can easily see, the dynamics of Ut 
for o < o f and also at o = <7f remains the same as discussed in the earlier section. At each fixed point there will be 
an equilibrium elongation 0(0) and a corresponding stress S = U* k5(<t) develops in the system (bundle). This 6(a) 
can be easily expressed in terms of U*(a). This requires the evaluation of a*, the internal stress per fiber developed 
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FIG. 6: Schematic stress (S')-strain (8) curve of the bundle (shown by the solid line), following Eq. I19H . with the fiber strength 
distribution llri^ (as shown in Fig. 5). Note, the model gives analytic solution for the full nonlinear stress-strain relationship 
of the bundle (disordered solid), including its failure (fracture) stress or strain. 



at the fixed point, corresponding to the initial (external) stress a (= F/N) per fiber applied on the bundle when all 
the fibers were intact. Expressing the effective stress a* per fiber in terms of U*(a), one can write from i|17|) 



U*(a) 

for a > ol- Consequently, 



<t L 1 - ff* 



k5{<j) = a* = 1 - (1 - cr L )U*(a). 

It may be noted that the internal stress a* (= a c /U*(a c ) is universally equal to 1/2 (independent of ox; from (|18(l ) 
at the failure point a = o c of the bundle. This finally gives the stress-strain relation for the model : 

{k5, < a < <jl 

k6(1-k5)/(1-<t l ), <T L <a<a f . (19) 
0, a > <T f 

This stress-strain relation is schematically shown in Fig. 6, where the initial linear region has slope k (the force 
constant of each fiber). This Hooke's region for stress S continues up to the strain value S — ctl/k, until which no 
fibers break (U*(a) = 1). After this, nonlinearity appears due to the failure of a few of the fibers and the consequent 
decrease of U*(<r) (from unity). It finally drops to zero discontinuously by an amount o*[/*(o/) = 1/ [4(1 — (Tl)\ = erf 
at the breaking point a = <jf or 8 = <j*Jk = 1/2/t for the bundle. This indicates that the stress drop at the final 
failure point of the bundle is related to the extent (o\l) of the linear region of the stress-strain curve of the same 
bundle. 



3. Strength of The Local Load Sharing Fiber Bundles 

So far, we studied models with fibers sharing the external load equally. This type of model shows (both analytically 
and numerically) existence of a critical strength (non zero aj) of the macroscopic bundle 01 beyond which it collapses. 
The other extreme model, i.e., the local load sharing model has been proved to be difficult to tackle analytically. 

It is clear, however, that the extreme statistics comes into play for such load load sharing models, for which the 
strength <jf — > as the bundle size (TV) approaches infinity. Essentially, for any finite load (a), depending on the fiber 
strength distribution, the size of the defect cluster can be estimated using Lifshitz argument (sec section 2.2) as In AT, 
giving the failure strength aj ~ l/(lniV) a , where the exponent a assumes a value appropriate for the model (see e.g., 
jjj)- If a fraction / of the load of the failed fiber goes for global redistribution and the rest (fraction 1 — /) goes to 
the fibers neighboring to the failed one, then we see (see Pradhan et al [T(|) that there is a crossover from extreme to 
self-averaging statistics at a finite value of /. 
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FIG. 7: A typical fracture surface: it has the self-afflne scaling property z(Xx, Xy) — X"z(x,y) where the roughness exponent 
£ has some universal value (e.g., £ ~ 0.8 for (2 + l)-dimensional fractured surface). 



III. TWO FRACTAL OVERLAP MODEL OF EARTHQUAKE AND ITS STATISTICS 

Overlapping fractals form a whole class of models to simulate earthquake dynamics. These models are motivated 
by the observation that a fault surface, like a fractured surface, is a fractal object Consequently a fault may 

be viewed as apair of overlapping fractals. Fractional Brownian profiles have been commonly used as models of 
fault surfaces [3|- In that case the dynamics of a fault is represented b y o ne Brownian profile drifting on another 
and each intersection of the two profiles corresponds to an earthquake [20|. However the simplest possible model 
of a fault — from the fractal point of view — was proposed by Chakrabarti and Stinchcombe |17|. This model is 
a schematic representation of a fault by a pair of dynamically overlapping Cantor sets. It is not realistic but, as a 
system of overlapping fractals, it has the essential feature. Since the Cantor set is a fractal with a simple construction 
procedure, it allows us to study in detail the statistics of the overlap of one fractal object on another. The two fractal 
overlap magnitude changes in time as one fractal moves over the other. The overlap (magnitude) time series can 
therefore be studied as a model time series of earthquake avalanche dynamics |l5j . 

The statistics of overlaps between two fractals is not studied much yet, though their knowledge is often required in 
various physical contexts. It has been established recently that since the fractured surfaces have got well-characterized 
sclf-affinc properties, the distribution of the elastic energies released during the slips between two fractal surfaces 
(earthquake events) may follow the overlap distribution of two self-similar fractal surfaces 0, 0] . Chakrabarti and 
Stinchcombe |17| had shown analytically by renormalization group calculations that for regular fractal overlap (Cantor 
sets and carpets) the contact area distribution p(s) follows a simple power law decay: 

p(s) ~ s-T; 7 = 1. (20) 

In this so called Chakrabarti-Stinchcombe model 0, the solid-solid contact surfaces of both the earth's crust and 
the tectonic plate are considered as average self-similar fractal surfaces. We then consider the distribution of contact 
areas, as one fractal surface slides over the other. We relate the total contact area between the two surfaces to be 
proportional to the elastic strain energy that can be grown during the sticking period, as the solid-solid friction force 
arises from the elastic strains at the contacts between the asperities. We then consider this energy to be released as 
one surface slips over the other and sticks again to the next contact or overlap between the rough surfaces. Since the 
two fractured surfaces are well known fractals, with established (self-affine) scaling properties (see Fig. Considering 
that such slips occur at intervals proportional to the length corresponding to that area, we obtain a power law for the 
frequency distribution of the energy releases. This compares quite well with the Gutenberg-Richter law. 

In order to proceed with the estimate of the number density n(e) of earthquakes releasing energy e in our model, 
we first find out the distribution p(s) of the overlap or contact area s between two self-similar fractal surfaces. We 
then relate s with e and the frequency of slips as a function of s, giving finally n(e). To start with a simple problem of 
contact area distribution between two fractals, we first take two Cantor sets [Tj] to model the contact area variations 
of two (nonrandom and self-similar) surfaces as one surface slides over the other. Figure [HJa) depicts structure in 
such surfaces at a scale which corresponds to only the second generation of iterative construction of two displaced 



12 




FIG. 8: Schematic representations of a portion of the rough surfaces of the earth's crust and the supporting (moving) tectonic 
plate, (b) The one dimensional projection of the surfaces form Cantor sets of varying contacts or overlaps as one surface slides 
over the other. 




La's tfc) ii-: 

FIG. 9: (a) Two cantor sets (in their first generation) along the axes r and r — r'. (b) This gives the overlap si(r) along the 
diagonal, (c) The corresponding density pi(s) of the overlap s at this generation. 

Cantor sets, shown in Fig. Elb). It is obvious that with successive iterations, these surfaces will acquire self-similarity 
at every length scale, when the generation number goes to infinity. We intend to study the distribution of the total 
overlap s (shown by the shaded regions in Fig. Elb)) between the two Cantor sets, in the infinite generation limit. 

Let the sequence of generators Gi define our Cantor sets within the interval [0,1]: Go = [0,1], G\ = RGo = 
[0, a] (J [6, 1] (i.e., the union of the intervals [0, a] and [6,1]), ■■■ , Gi+i — RGi, ... . If we represent the mass density 
of the set Gi by Di(r), then Di(r) = 1 if r is in any of the occupied intervals of Gi, and Di(r) = elsewhere. 
The required overlap magnitude between the sets at any generation I is then given by the convolution form s;(r) = 
/ dr'Di(r')Di(r — r'). This form applies to symmetric fractals (with Di(r) = Di(—r)); in general the argument of the 
second D; should be Di(r + r'). 

One can express the overlap integral s% in the first generation by the projection of the shaded regions along the 
vertical diagonal in Fig. 8(a). That gives the form shown in Fig. 8(b). For a = b < i, the nonvanishing Si(r) regions 
do not overlap, and are symmetric on both sides with the slope of the middle curve being exactly double those on the 
sides. One can then easily check that the distribution pi(s) of overlap s at this generation is given by Fig. 8(c), with 
both c and d greater than unity, maintaining the normalisation of the probability pi with cd — 5/3. The successive 
generations of the density pi(s) may thefore be represented by Fig. |3 where 

Pl+1 (s) = R Pl {s) = - Pl (-J + —pi I- . (21) 

In the infinite generation limit of the renormalisation group (RG) equation, if p*(s) denotes the fixed point distribution 




(a) $ m s (c) s (A) $ 



FIG. 10: The overlap densities p(s) at various generations of the Cantor sets at the zeroth (a), first (b), second (c) and at the 
infinite (or fixed point) (d) generations. 
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FIG. 11: For two (regular) Cantor sets, one moving uniformly over the other (with periodic boundary condition), the total 
measure of the (shaded region in Fig. |Ub) contributing to the) overlap s has the time variation as shown here (for n = 4). 

such that p*(s) = Rp*(s), then assuming p*(s) ~ s~ 7 /5(s), one gets (d/5)c 7 + (4d/5)(c/2) 7 = 1. Here p(s) represents 
an arbitrary modular function, which also includes a logarithmic correction for large s. This agrees with the above 
mentioned normalisation condition cd — 5/3 for the choice 7 = 1. This result for the overlap distribution l|2U|) 



p*(s) = p(s) 



7=1, 



is the general result for all cases that we have investigated and solved by the functional rescaling technique (with the 
logs correction for large s, renormalising the total integrated distribution). 

The above study is for continuous relative motion of one Cantor set over the other. Study of the time (t) variation 
of contact area (overlap) s(t) between two well-characterized fractals having the same fractal dimension as one fractal 
moves over the other with constant velocity, with discrete (minimum element in that generation) steps has been 
studied, for finite generations [18|. Bhattacharyya [l8j studied this overlap distribution for two Cantor sets with 
periodic boundary conditions and each having dimension log 2 / log 3 (see Fig. Illf) . It was shown, using exact counting, 
that if s = 2 n ~ k (n is the generation number) then the probability p(s) to get an overlap s is given by a binomial 
distribution 



~(<~>n — k \ 



p(2 



n 

n — k 



exp(— r 2 /n); r — » 0, 



(22) 



where r 2 = [| — &)] 2 - Expressing therefore r by logs near the maxima of p(s), one can again rewrite l12L'l) as 



p(s) ~ exp - 



(logs)' 



(23) 



Noting that p(s)d(logs) ~ p(s)ds, we find p(s) ~ s~ 7 , j = 1, as in H2U|) . as the binomial or Gaussian part becomes 
a very weak function of s as n — * oo [2l|. It may be noted that this exponent value 7 = 1 is independent of the 
dimension of the Cantor sets considered (here log 2/ log 3) or for that matter, independent of the fractals employed. 
It also denotes the general validity of (|2*U|> even for disordered fractals, as observed numerically @. 

Identifying the contact area or overlap s between the self-similar (fractal) crust and tectonic plate surfaces as the 
stored elastic energy E released during the slip, the distribution (|20(l . of which a derivation is partly indicated here, 
reduces to the Gutenberg-Richter law © observed. 



IV. SUMMARY AND DISCUSSIONS 



Unlike the elastic constants (e.g., Y in eqn. @J) of a solid, which can be estimated from the interatomic interactions 
and lattice structures of a solid and while the effect of disorder on them can be accommodated using simple analytic 
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formulas, the fracture strength (ay) of a solid cannot be estimated easily from there and such estimates are orders 
of magnitude higher than those observed. The reason is that cracks nucleate around a defect or disorder in the solid 
and the variation of ay with the defect size I can be precisely estimated in a brittle solid using the Griffith's formula 
@ in Sec. 2.1. 

For disordered solids, the failure strength distribution F(a) depends on the sample volume and is given by the 
extreme statistics: Gumbel0] or Weibull Q type (see Fig. 2 for their generic behavior). The average strength of 
a (finite volume) sample can be estimated from a modified Griffith-like formula Q. When percolation correlation 
length exceeds the Lifshitz length, for large disorder, self-averaging statistics takes over (see Sec. 2.2) and the average 
strength, given by (J7|), becomes precisely defined (even for infinite system size) and a/ becomes volume independent 
(as in the equal load sharing fiber bundle model). 

The inherent mean-field nature of the equal load sharing models (discussion in Sec. 2.3) enables us to construct 
recursion relations (Eq. © for example) which captures essentially all the intriguing features of the failure dynamics. 
Though we have identified O = U*(a) — U*(<jf) cx (a/ — a) 13 as the order parameter (with exponent (3 — 1/2) for 
the continuous transition in such models, unlike in the conventional phase transitions it does not have a real- valued 
existence for a > <jf. The 'type' of phase transition in such models has been a controversial issue. Earlier it was 
suggested to be a first order phase transition, because the the surviving fraction of fibers has a discontinuity at the 
breakdown point of the bundles. However, as the susceptibility shows divergence (x oc (ay — a)~ 7 ;7 = 1/2) at the 
breakdown point, the transition has been later identified to be of second order |ldj| . The dynamic critical behavior of 
the these models and the universality of the exponent values are straightforward. Here, divergence of relaxation time 
(r) at the critical point (r cx (ay — a)~ Q ; a — 1/2) indicates 'critical slowing' of the dynamics which is characteristic 
of conventional critical phenomena. At the critical point, one observes power law decay of the surviving fraction in 
time (Ut(<Jf) oc t~ s ;5 = 1). We demonstrated the universality of the failure behavior near a = a/, for three different 
distributions: uniform, linearly increasing and linearly decreasing distributions of fiber strength. The critical strengths 
of the bundles differ in each case: a/ = 1/4, and 4/27 respectively for these three distributions. However, the 

critical behavior of the order parameter O, susceptibility x, relaxation time r and of the time decay at ay, as given 
by the exponents /3, 7, a and S remain unchanged: a=l/2 = /3 = 7 and 6 = 1 for all three distributions. 

The model also shows realistic nonlinear deformation behavior with a shifted (by ul , away from the origin) uniform 
distribution of fiber strengths (see Sec. 2.3.2). The stress-strain curve for the model clearly shows three different 
regions: elastic or linear part (Hooke's region) when none of the fibers break (U*(a) = 1), plastic or nonlinear part due 
to the successive failure of the fibers (C/*(a) < 1) and then finally the stress drops suddenly (due to the discontinuous 
drop in the fraction of surviving fibers from U*(af) to zero) at the failure point a/ = 1/ [4(1 — ax)]. This nonlinearity 
in the response (stress-strain curve in Fig. 6) results from the linear response of the surviving fibers who share the 
extra load uniformly. The local load sharing bundles (see Sec. 2.3.3) on the other hand show 'zero' critical strength 
as the bundle size goes to infinity in one dimension (extreme statistics takes over). It is not clear at this stage if, in 
higher dimensions, LLS bundles are going to have non-zero critical strength. In any case, the associated dynamics of 
failure of these higher dimensional bundles with variable range load transfer should be interesting. 

We believe, the elegance and simplicity of the model, its common-sense appeal, the exact solubility of its critical 
behavior in the mean field (ELS) limit, its demonstrated universality, etc, would promote the model eventually to a 
level competing with the Ising model of magnetic critical behavior. 

As emphasized already, we consider the physicist's identification of the Gutenberg- Richter law (3 ) as an extremely 
significant one in geophysics. Like the previous attempts [l4l Il5l [l q . the model developed here |17| captures this 
important feature in its resulting statistics. Here, the established self-similarity of the fault planes are captured using 
fractals, Cantor sets in particular. Hence we consider, in Sec. 3, this 'Chakrabarti-Stinchcombe' model 0]> where 
one Cantor set moves uniformly over another similar set (with periodic boundary conditions). The resulting overlap 
s (meaning the set of real numbers common in both the Cantor sets) changes with time: see, for example, Fig. II II for 
a typical time variation overlap s for n — 4. The number density of such overlaps seem to follow a Gutenberg-Richter 
type law 1|20|) . Judging from the comparisons of the exponent values a in and 7 in (|20|l . the model succeeds 
at least as well as the earlier ones. More importantly, our model incorporates both the geologically observed facts: 
fractal nature of the contact surfaces of the crust and of the tectonic plate, and the stick-slip motion between them. 
However, the origin of the power law in the quake statistics here is the self-similarity of the fractal surfaces, and 
not any self-organisation directly in their dynamics. In fact, the extreme non-linearity in the nature of the crack 
propagation is responsible for the fractal nature of the rough crack surfaces of the crust and the tectonic plate. This 
in turn leads here to the Gutenberg-Richter like power law in the earthquake statistics. 

Acknowledgments: The author thanks M. Acharyya, K. K. Bardhan, L. G. Bcnguigui, P. Bhattacharyya, A. 
Chatterjee, D. Chowdhury, M. K. Dey, A. Hansen, S. S. Manna, S. Pradhan, P. Ray, D. Stauffer and R. B. Stinchcombe 



15 



for collaborations at different stages. 



[1] B. Lawn, Fracture of Brittle Solids, Cambridge Univ. Press, Cambridge (1993). 

[2] B. K. Chakrabarti and L. G. Benguigui, Statistical Physics of Fracture and Breakdown in Disorder Systems, Oxford Univ. 
Press, Oxford (1997). 

[3] H.J. Herrmann and S. Roux (Eds.) Statistical Models for the Fracture of Disordered Media, Elsevier, Amsterdam (1990); 

M. Sahimi, Heterogeneous Materials, Vol. II, Springer, New York (2003). 
[4] D. Stauffer and A. Aharony, Introduction to Percolation Theory, Taylor and Francis, London (1992). 
[5] A. A. Griffith, Phil. Trans. Roy. Soc. Lond. A 221 163 (1920). 
[6] F. T. Pierce, J. Textile Inst. 17, T355-368 (1926). 

[7] H. E. Daniels, Proc. R. Soc. London A 183 405 (1945); S. L. Phoenix, SIAM J. Appl. Math. 34 227 (1978); Adv. Appl. 

Prob. 11 153 (1979). 
[8] S. Pradhan and B. K. Chakrabarti, Int. J. Mod. Phys. B 17 5565 (2003). 

[9] P. C. Hemmer and A. Hansen, J. Appl. Mech. 59 909 (1992); M. Kloster, A. Hansen and P. C. Hemmer, Phys. Rev. E 56 

2615 (1997); S. Pradhan, A. Hansen and P. C. Hemmer, Phys. Rev. Lett. 95 125501 (2005); F. Raischel, F. Kun and H. 

J. Herrmann, cond-mat/0601290 (2006). 
[10] S. Pradhan and B. K. Chakrabarti, Phys. Rev. E 65 016113 (2001); S. Pradhan, P. Bhattacharyya and B. K. Chakrabarti, 

Phys. Rev. E 66 016116 (2002); P. Bhattacharyya, S. Pradhan and B. K. Chakrabarti, Phys. Rev. E 67 046122 (2003). 
[11] R. C. Hidalgo, F. Kun and H. J. Herrmann, Phys. Rev. E 64 066122 (2001); S. Pradhan, B. K. Chakrabarti and A. Hansen, 

Phys. Rev. E 71 036149 (2005). 
[12] B. Gutenberg and C. F. Richter, Seismicity of the Earth and Associated Phenomena, Princeton Univ. Press, Princeton, 

N.J. (1954). 

[13] L. Knopoff, Proc. Natl. Acad. Sci. USA 97 11880 (2000); Y. Y. Kagan, Physica D 77 160 (1994); C. H. Scholz, The 
Mechanics of Earthquake and Faulting, Cambridge Univ. Press, Cambridge (1990); B. V. Kostrov and S. Das, Principles 
of Earthquake Source Mechanics, Cambridge Univ. Press, Cambridge (1988). 

[14] R. Burridge and L. Knopoff, Bull. Seis. Soc. Am. 57 341-371 (1967). 

[15] J. M. Carlson and J. S. Langer, Phys. Rev. Lett. 62 2632-2635 (1989); J. M. Carlson, J. S. Langer and B. E. Shaw, 
Rev. Mod. Phys. 66 657-670 (1994); G. Ananthakrishna and H. Ramachandran in Nonlmearity and Breakdown in Soft 
Condensed Matter, Eds. K. K. Bardhan, B. K. Chakrabarti and A. Hansen, LNP 437, Springer Verlag, Heidelberg (1994); 
T. Mori and H. Kawamura, Phys. Rev. Letts., 94 058501 (2005). 

[16] P. Bak, How Nature works, Oxford Univ. Press, Oxford (1997). 

[17] B. K. Chakrabarti and R. B. Stinchcombe, Physica A 270 27 (1999); S. Pradhan, B. K. Chakrabarti, P. Ray and M. K. 

Dey, Phys. Scripta T106 77 (2003). 
[18] P. Bhattacharyya, Physica A 348 199 (2005). 

[19] D. J. Bergman and D. Stroud, in Solid State Physics, 46 Eds. H. Ehrenreich and D. Turnbull, Academic Press, New York, 
p.147 (1992). 

[20] V. de Rubeis, R. Hallgass, V. Loreto, G. Paladin, L. Pietronero and P. Tosi, Phys. Rev. Lett. 76 2599 (1996). 

[21] B. K. Chakrabarti and A. Chatterjee, in Proc. The Seventh International Conference on Vibration Problems ICOVP-2005, 

Istanbul, Ed. E. Inan (Springer, 2006), arXiv:cond-mat/0512136 P. Bhattacharyya, A. Chatterjee and B. K. Chakrabarti, 

arXiv:physics/0510038 



